home *** CD-ROM | disk | FTP | other *** search
/ SGI Hot Mix 17 / Hot Mix 17.iso / HM17_SGI / research / lib / factorial.pro < prev    next >
Text File  |  1997-07-08  |  2KB  |  67 lines
  1. ;$Id: factorial.pro,v 1.3 1997/01/15 03:11:50 ali Exp $
  2. ;
  3. ; Copyright (c) 1994-1997, Research Systems, Inc.  All rights reserved.
  4. ;       Unauthorized reproduction prohibited.
  5. ;+
  6. ; NAME:
  7. ;       FACTORIAL
  8. ;
  9. ; PURPOSE:
  10. ;       This function computes the factorial N! as the double-precision 
  11. ;       product, (N) * (N-1) * (N-2) * ...... * 3 * 2 * 1. 
  12. ;
  13. ; CATEGORY:
  14. ;       Special Functions.
  15. ;
  16. ; CALLING SEQUENCE:
  17. ;       Result = Factorial(n)
  18. ;
  19. ; INPUTS:
  20. ;       N:    A non-negative scalar of type integer, float or double.
  21. ;
  22. ; KEYWORD PARAMETERS:
  23. ;       STIRLING:    If set to a non-zero value, Stirling's asymptotic 
  24. ;                    formula is used to approximate N!. 
  25. ;
  26. ; EXAMPLE:
  27. ;       Compute 20! with and without Stirling's asymptotic formula.
  28. ;         result_1 = factorial(20, /stirling)
  29. ;         result_2 = factorial(20)
  30. ;       
  31. ;       Result_1 and result_2 should be 2.4227869e+18 and 2.4329020e+18
  32. ;       respectively.
  33. ;
  34. ; REFERENCE:
  35. ;       ADVANCED ENGINEERING MATHEMATICS (seventh edition)
  36. ;       Erwin Kreyszig
  37. ;       ISBN 0-471-55380-8
  38. ;
  39. ; MODIFICATION HISTORY:
  40. ;       Written by:  GGS, RSI, November 1994
  41. ;-
  42.  
  43. function factorial, n, stirling = stirling
  44.  
  45.   ;Computes N! as (N) * (N-1) * (N-2) * ...... * 3 * 2 * 1
  46.   ;Test example: 20! = 2.4329020e+18
  47.   ;Use NR_MACHAR(/DOUBLE) to determine largest floating point number.
  48.   ;Stirling's formula: N! = SQRT(2.0d*!PI*N) * (N/EXP(1.0d))^(N+0.0d)
  49.  
  50.   on_error, 2
  51.  
  52.   if n lt 0 then $
  53.     message, 'n must be a non-negative scalar.'
  54.  
  55.   ;Computes N! as (N) * (N-1) * (N-2) * ...... * 3 * 2 * 1
  56.   if keyword_set(stirling) eq 0 then begin 
  57.     fact = 1.0d
  58.     for k = n+0.0d, 1d, -1d do begin
  59.       fact = fact * k
  60.     endfor
  61.   endif else $ ;Approximate N! using Stirling's formula.
  62.     fact = sqrt(2.0d * !pi * n) * (n / exp(1.0d))^(n+0.0d)
  63.  
  64.   return, fact
  65.  
  66. end
  67.